structure_analysis
Amanda Stahlke
2/28/2020
Structure Analysis
Examined "global" ancestry for all available samples. Old replicates and some suspected contamination/switching have been removed. The masterpopmap has ALL of these samples.
Admixture model, LD=0, allele freqs are uncorrelated, alpha is inferred for each popualtion.
Read in Results and set up popmap/poporder for visualization and downstream analysis.
sfiles <- list.files(path="Results/",
pattern=glob2rx("outfile_*f"),full.names=TRUE)
slist <- readQ(sfiles,filetype="structure",
indlabfromfile = TRUE,
readci = TRUE)
popmap <- left_join(data.frame(sampID = rownames(slist[[1]])),
masterpopmap)## Joining, by = "sampID"popID <- data.frame(popID = popmap$popID)
popID <- data.frame(lapply(popID, as.character), stringsAsFactors=FALSE)
pop_order <- sample_meta.dat[,c("popID", "Order")]
pop_order <- arrange(pop_order, Order)[,1]Sample size for just Structure analysis
samplesize <- left_join(popmap, sample_meta.dat) %>%
count(popID) %>%
merge(sample_meta.dat)## Joining, by = "popID"samplesize## popID n
## 1 11CO 12
## 2 12AZ 15
## 3 15TX 16
## 4 16TX 18
## 5 18TX 17
## 6 19TX 15
## 7 1WY 17
## 8 20NM 17
## 9 21NM 14
## 10 26NM 18
## 11 27NM 20
## 12 28NM 17
## 13 2CO 19
## 14 31NV 7
## 15 32UT 20
## 16 33NV 9
## 17 34UT 16
## 18 37CR 5
## 19 38CR 11
## 20 39CR 12
## 21 41CR 16
## 22 43GR 15
## 23 44GR 16
## 24 46CH 24
## 25 47TX 25
## 26 49TX 27
## 27 4CO 11
## 28 50TX 8
## 29 51TX 14
## 30 52TX 21
## 31 5CO 19
## 32 6KS 10
## 33 8OK 20
## 34 CARINA_LAB 14
## 35 SUB_LAB 17
## Name Latitude
## 1 Wilkinson 37.34
## 2 Big Bend 35.12
## 3 NE Post\n 33.32
## 4 Lake JB Thomas 32.61
## 5 Orla 31.49
## 6 Aspermont 33.17
## 7 Lovell 44.86
## 8 Malaga\n 32.22
## 9 Artesia Wildlife Reserve 32.98
## 10 Lake Sumner 34.15
## 11 Roswell E 33.40
## 12 Tucumcari Lake 35.19
## 13 Fountain Creek 38.34
## 14 Virgin River - Gold Butte 36.69
## 15 St. George 37.07
## 16 Lake Mead, Stewarts Pt 36.38
## 17 Delta 39.23
## 18 Rethimno, Crete, Greece 35.37
## 19 Crete, GR 35.42
## 20 Panaramnos, Crete, Greece 35.42
## 21 Plakias, Crete, Greece 35.19
## 22 Posidi, Greece 39.97
## 23 Delta Aksiou, Greece 40.55
## 24 Bitun, China 47.30
## 25 Rio Grande Village - from T. cheninsis x T. ramosissima hybrids 29.18
## 26 Santa Elana 29.16
## 27 Adobe Reservoir (Blue Lake)\n 38.26
## 28 Presidio Hwy 29.34
## 29 CO Canyon Boat Ramp 29.34
## 30 Tornillo 31.44
## 31 SE corner 37.70
## 32 W Finney County 37.99
## 33 Guymon 36.70
## 34 D. carinata Lab Culture\n NA
## 35 D. sublineata Lab Culture\n NA
## Longitude Date.Collected Order Grp cat
## 1 -104.16 9/18/14 7 carinulata expansion
## 2 -114.64 9/18/14 11 carinulata expansion
## 3 -101.26 9/26/14 27 potential_hybrid_zone expansion
## 4 -101.22 9/26/14 30 potential_hybrid_zone original_release
## 5 -103.48 9/27/14 32 potential_hybrid_zone original_release
## 6 -100.24 9/27/14 28 potential_hybrid_zone expansion
## 7 -108.18 9/3/14 2 carinulata original_release
## 8 -104.08 9/27/14 31 potential_hybrid_zone original_release
## 9 -104.44 10/7/14 29 potential_hybrid_zone original_release
## 10 -104.48 9/29/14 25 potential_hybrid_zone original_release
## 11 -104.41 9/29/14 26 potential_hybrid_zone original_release
## 12 -103.69 10/6/14 24 potential_hybrid_zone original_release
## 13 -104.61 9/2/14 4 carinulata original_release
## 14 -114.26 10/15/14 9 carinulata expansion
## 15 -113.58 10/15/14 8 carinulata expansion
## 16 -114.40 10/15/14 10 carinulata expansion
## 17 -112.93 10/15/14 3 carinulata original_release
## 18 24.47 7/2/15 18 native_elongata native
## 19 24.69 7/4/15 17 native_elongata native
## 20 24.68 7/4/15 16 native_elongata native
## 21 24.40 7/4/15 20 native_elongata native
## 22 23.37 7/6/15 15 native_elongata native
## 23 22.74 7/6/15 14 native_elongata native
## 24 87.75 7/3/16 1 carinulata native
## 25 -102.96 8/3/17 36 potential_hybrid_zone original_release
## 26 -103.60 8/3/17 37 potential_hybrid_zone original_release
## 27 -103.25 9/4/14 5 carinulata expansion
## 28 -104.07 8/3/17 34 potential_hybrid_zone original_release
## 29 -104.06 8/3/17 35 potential_hybrid_zone original_release
## 30 -106.09 8/3/17 33 potential_hybrid_zone expansion
## 31 -103.42 8/21/14 6 carinulata expansion
## 32 -101.08 8/8/14 22 potential_hybrid_zone expansion
## 33 -101.55 9/16/14 23 potential_hybrid_zone expansion
## 34 NA 12 lab lab
## 35 NA 13 lab labwrite.table(samplesize, "pop_metadata_samplesize.tsv", sep = "\t", row.names = F, quote = F)What is the best K?
sr1 <- summariseQ(tabulateQ(slist))
p <- evannoMethodStructure(data=sr1,exportplot=F,returnplot=T,returndata=F,basesize=12,linesize=0.7)
grid.arrange(p)
evannoMethodStructure(data=sr1,
exportplot=T,
exportpath = "plots",
outputfilename = "global_deltaK",
returndata=F,
height = 10, width = 9)## plots/global_deltaK.png exported.K=2 and K=3 are super stable. D. carinata appears as a mix between sublineata and elongata. I am surprised to see, for the first time in these samples, indications of admixture between carinlata and sublineata along the Rio Grande.
As usual, K=4 reflects lab populations, elongata in the native range, and hints at carinulata substructure.
k4 <- plotQ(alignK(slist[c(41:50)]),
imgoutput="join",
returnplot=T,exportplot=F,basesize=11,
clustercol=cbbPalette,
grplab = popID,
ordergrp = T,
grplabangle=-90)
grid.arrange(k4$plot[[1]])
Code to save a nice version of all reps of K=4 for the manuscript.
fn1 <- function(x) attr(x,"k")
spnames <- rep("", 10)
plotQ(alignK(slist[c(41:50)]),
imgoutput="join",
exportpath = "plots",
outputfilename = "Joined-K4",
exportplot=T,
showyaxis=F,
showsp = T,
splab=rep("", 10),
clustercol=cbbPalette,
subsetgrp = pop_order,
height = .005,
ordergrp = T,
grplab = popID,
linepos = 1,
divtype = 1,
grplabangle=90,
grplabpos = .8,
grplabspacer = 0,
basesize=14,
grplabsize = 15,
width = 150,
grplabheight=150)## Drawing plot ...## plots/Joined-K4.png exported.Assign species identity to each cluster
Extract one q file that reflects known species identity
global_k4_q <- slist[[41]]
global_k4_q$sampID <- row.names(global_k4_q)
## which cluster is carinulata?
gr_samples <- masterpopmap[grep( "CR", masterpopmap[,2]), ]
gr_q <- semi_join(global_k4_q, gr_samples)
max_cl_gr <- sapply(gr_q[,1:4], max)
elongata_cluster <- which(colnames(gr_q) == colnames(gr_q[which(max_cl_gr==1)]))
elongata_cluster## [1] 1## which cluster is elongata?
ch_samples <- masterpopmap[grep( "CH", masterpopmap[,2]), ]
ch_q <- semi_join(global_k4_q, ch_samples)
max_cl_ch <- sapply(ch_q[,1:4], max)
carinulata_cluster <- which(colnames(ch_q) == colnames(ch_q[which(max_cl_ch==1)]))
carinulata_cluster## [1] 3## which cluster is sublineata?
sub_samples <- masterpopmap[grep( "SUB_LAB", masterpopmap[,2]), ]
sub_q <- semi_join(global_k4_q, sub_samples)
max_cl_sub <- sapply(sub_q[,1:4], max)
sublineata_cluster <- which(colnames(sub_q) == colnames(sub_q[which(max_cl_sub==1)]))
## which cluster is carinata?
carina_samples <- masterpopmap[grep( "CARINA_LAB", masterpopmap[,2]),]
carina_q <- semi_join(global_k4_q, carina_samples)
max_cl_cara <- sapply(carina_q[,1:4], max)
carinata_cluster <- which(colnames(carina_q) == colnames(carina_q[which(max_cl_cara==1)]))
cls <- c(carinata_cluster, carinulata_cluster, sublineata_cluster, elongata_cluster)
spp <- c("carinata", "carinulata", "sublineata", "elongata")
clusters <- as.data.frame(cbind(spp = spp[order(cls)], cluster = seq(1,4)))
colnames(global_k4_q)[1:4] <- as.character(clusters$spp)
q_popmap <- merge(global_k4_q, popmap)
head(q_popmap)## sampID elongata sublineata carinulata carinata popID
## 1 D1 1.000 0.000 0.00 0.000 38CR
## 2 D12 1.000 0.000 0.00 0.000 41CR
## 3 D13 1.000 0.000 0.00 0.000 41CR
## 4 D15_rep 0.970 0.000 0.03 0.000 41CR
## 5 D16 1.000 0.000 0.00 0.000 41CR
## 6 D17 0.998 0.001 0.00 0.001 43GRlegend <- merge(legend_colors, clusters)
clustercols <- arrange(legend, cluster) %>%
select(cbbPalette)
clustercols <- clustercols$cbbPalette
legendlab <- arrange(legend, cluster) %>%
mutate(label = paste("D. ", spp))
legendlab <- legendlab$label[1:4]Code to save a nice version of one rep of K=4 for the manuscript.
plotQ(slist[41],
exportpath = "plots",
outputfilename = "MainTextFigure_K4",
exportplot=T,
showyaxis=T,
showsp = F,
clustercol=clustercols,
subsetgrp = pop_order,
sortind="all",
grplab = popID,
linepos = 1,
linesize=0.8,
pointsize = 5,
divtype = 1,
divsize=1.1,
grplabangle= 90,
grplabpos = .9,
grplabspacer = 0,
grplabheight=35,
grplabsize = 10,
height = .05,
width=100,
showlegend = T, legendkeysize=30,legendtextsize=30, legendpos="left", legendlab=legendlab,
# titlelab="Global Population Structure K=4",
# titlesize = 20,
basesize=30)Code to save a nice version of K=2-4 for the manuscript.
# custom strip panel label showing k only
fn1 <- function(x) attr(x,"k")
spnames <- paste0("K=",sapply(slist,fn1))
# Combine K=2, K=3, K=4
k2k3k4 <- c(21,31,41)
plotQ(slist[k2k3k4],
imgoutput="join",
exportpath = "plots",
outputfilename = "Joined-K2-K3-K4",
exportplot=T,
showyaxis=F,
clustercol=clustercols,
subsetgrp = pop_order,
ordergrp = T,
grplab = popID,
linepos = 1,
divtype = 1,
# divsize = 1.2,
grplabangle= 90,
grplabpos = .8,
grplabspacer = 0,
basesize=14,
showsp = F,
grplabsize = 6,
width = 30,
grplabheight=40)Exploratory visualizations that are not presented in the manuscript
Multi-line K=4 to visulize individuals
plotQMultiline(slist[41],
showticks=T,
clustercol=cbbPalette,
grplab = popID,
ordergrp=TRUE,
spl = 25,
showindlab = T,
showyaxis = T,
barsize = .9,
subsetgrp = pop_order,
basesize=11,
exportpath = "plots",
outputfilename = "multilineK4",
useindlab = T)## Drawing plot ...## plots/multilineK4-1.png exported.## Drawing plot ...## plots/multilineK4-2.png exported.## Drawing plot ...## plots/multilineK4-3.png exported.All reps of each run to visualize results across K values
k2 <- plotQ(alignK(slist[c(21:30)]),
imgoutput="join",
returnplot=T,exportplot=F,basesize=11,
clustercol=cbbPalette,
subsetgrp = pop_order,
ordergrp = T,
grplab = popID,
linepos = 1,
grplabangle=-45,
grplabpos = 0,
grplabheight=20)
grid.arrange(k2$plot[[1]])
k3 <- plotQ(alignK(slist[c(31:40)]),
imgoutput="join",
returnplot=T,exportplot=F,basesize=11,
clustercol=cbbPalette,
ordergrp = T,
grplab = popID,
grplabangle=-90)
grid.arrange(k3$plot[[1]])
k5 <- plotQ(alignK(slist[c(51:60)]),
imgoutput="join",
returnplot=T,exportplot=F,basesize=11,
clustercol=cbbPalette,
grplab = popID,
ordergrp = T,
grplabangle=-90)
grid.arrange(k5$plot[[1]])
k6 <- plotQ(alignK(slist[c(61:70)]),
imgoutput="join",
returnplot=T,exportplot=F,basesize=11,
clustercol=cbbPalette,
grplab = popID,
ordergrp = T,
grplabangle=-90)
grid.arrange(k6$plot[[1]])
k7 <- plotQ(alignK(slist[c(71:80)]),
imgoutput="join",
returnplot=T,exportplot=F,basesize=11,
clustercol=cbbPalette,
grplab = popID,
grplabangle=-90)
grid.arrange(k7$plot[[1]])
k8 <- plotQ(alignK(slist[c(81:90)]),
imgoutput="join",
returnplot=T,exportplot=F,basesize=11,
# clustercol=cbbPalette,
grplab = popID,
grplabangle=-90)
grid.arrange(k8$plot[[1]])
Among individuals, how many pure? How many hybrids, and of what type?
Extract confidence intervals.
ci <- data.frame(attributes(slist[[41]])$ci)
colnames(ci) <- paste0(rep(colnames(global_k4_q)[1:4], each=2),
"_", rep(c("L", "H"), times = 4))
ci <- cbind(ci, global_k4_q)
ci <- merge(ci, masterpopmap)
head(ci)## sampID elongata_L elongata_H sublineata_L sublineata_H carinulata_L
## 1 D1 0.999 1.000 0 0.001 0.000
## 2 D12 0.998 1.000 0 0.001 0.000
## 3 D13 0.998 1.000 0 0.001 0.000
## 4 D15_rep 0.953 0.985 0 0.002 0.015
## 5 D16 0.999 1.000 0 0.001 0.000
## 6 D17 0.989 1.000 0 0.005 0.000
## carinulata_H carinata_L carinata_H elongata sublineata carinulata carinata
## 1 0.001 0 0.001 1.000 0.000 0.00 0.000
## 2 0.001 0 0.001 1.000 0.000 0.00 0.000
## 3 0.001 0 0.001 1.000 0.000 0.00 0.000
## 4 0.046 0 0.001 0.970 0.000 0.03 0.000
## 5 0.001 0 0.001 1.000 0.000 0.00 0.000
## 6 0.001 0 0.005 0.998 0.001 0.00 0.001
## popID
## 1 38CR
## 2 41CR
## 3 41CR
## 4 41CR
## 5 41CR
## 6 43GRwrite.csv(ci,"ancestry_confidenceintervals.csv",
quote = FALSE)Use "diagnostic" samples to define threshold.
diag_sublineata <- c("SUB_LAB", "47TX", "49TX", "50TX", "51TX", "52TX")
diagnostic <- data.frame(popID = diag_sublineata, spp=rep("sublineata", length(diag_sublineata)))
diag_elongata <- c("43GR", "44GR", "37CR", "38CR", "39CR", "41CR")
diagnostic <- rbind(diagnostic,
data.frame(popID = diag_elongata, spp=rep("elongata", length(diag_elongata))))
ci %>%
semi_join(diagnostic) %>%
arrange(desc(carinulata_L))## Joining, by = "popID"## sampID elongata_L elongata_H sublineata_L sublineata_H carinulata_L
## 1 N7 0.000 0.001 0.888 0.934 0.066
## 2 D24_rep 0.910 0.951 0.000 0.001 0.049
## 3 D15_rep 0.953 0.985 0.000 0.002 0.015
## 4 D1 0.999 1.000 0.000 0.001 0.000
## 5 D12 0.998 1.000 0.000 0.001 0.000
## 6 D13 0.998 1.000 0.000 0.001 0.000
## 7 D16 0.999 1.000 0.000 0.001 0.000
## 8 D17 0.989 1.000 0.000 0.005 0.000
## 9 D18 0.998 1.000 0.000 0.001 0.000
## 10 D19 0.998 1.000 0.000 0.001 0.000
## 11 D2 0.999 1.000 0.000 0.001 0.000
## 12 D20 0.998 1.000 0.000 0.001 0.000
## 13 D21 0.997 1.000 0.000 0.001 0.000
## 14 D22 0.998 1.000 0.000 0.001 0.000
## 15 D23 0.998 1.000 0.000 0.001 0.000
## 16 D27 0.998 1.000 0.000 0.001 0.000
## 17 D3 0.998 1.000 0.000 0.001 0.000
## 18 D33 0.998 1.000 0.000 0.001 0.000
## 19 D34 0.999 1.000 0.000 0.001 0.000
## 20 D35 0.998 1.000 0.000 0.001 0.000
## 21 D37 0.999 1.000 0.000 0.001 0.000
## 22 D38 0.998 1.000 0.000 0.001 0.000
## 23 D39 0.998 1.000 0.000 0.001 0.000
## 24 D4 0.998 1.000 0.000 0.001 0.000
## 25 D40 0.998 1.000 0.000 0.001 0.000
## 26 D41 0.998 1.000 0.000 0.001 0.000
## 27 D42 0.998 1.000 0.000 0.001 0.000
## 28 D43 0.998 1.000 0.000 0.001 0.000
## 29 D44 0.999 1.000 0.000 0.001 0.000
## 30 D46 0.997 1.000 0.000 0.001 0.000
## 31 D47 0.998 1.000 0.000 0.001 0.000
## 32 D48 0.999 1.000 0.000 0.001 0.000
## 33 D49 0.998 1.000 0.000 0.001 0.000
## 34 D5 0.999 1.000 0.000 0.001 0.000
## 35 D50 0.999 1.000 0.000 0.001 0.000
## 36 D51 0.998 1.000 0.000 0.001 0.000
## 37 D52 0.999 1.000 0.000 0.001 0.000
## 38 D53 0.999 1.000 0.000 0.001 0.000
## 39 D54 0.998 1.000 0.000 0.001 0.000
## 40 D55 0.998 1.000 0.000 0.001 0.000
## 41 D56 0.999 1.000 0.000 0.001 0.000
## 42 D57 0.998 1.000 0.000 0.001 0.000
## 43 D58 0.998 1.000 0.000 0.001 0.000
## 44 D59 0.998 1.000 0.000 0.001 0.000
## 45 D6 0.999 1.000 0.000 0.001 0.000
## 46 D60 0.998 1.000 0.000 0.001 0.000
## 47 D61 0.999 1.000 0.000 0.001 0.000
## 48 D62 0.998 1.000 0.000 0.001 0.000
## 49 D63 0.998 1.000 0.000 0.001 0.000
## 50 D64 0.998 1.000 0.000 0.001 0.000
## 51 D65_rep 0.970 1.000 0.000 0.012 0.000
## 52 D7 0.999 1.000 0.000 0.001 0.000
## 53 D71_rep 0.988 1.000 0.000 0.001 0.000
## 54 D72 0.997 1.000 0.000 0.001 0.000
## 55 D73 0.998 1.000 0.000 0.001 0.000
## 56 D75 0.998 1.000 0.000 0.001 0.000
## 57 D76 0.998 1.000 0.000 0.001 0.000
## 58 D77 0.998 1.000 0.000 0.001 0.000
## 59 D78 0.999 1.000 0.000 0.001 0.000
## 60 D79 0.998 1.000 0.000 0.001 0.000
## 61 D8 0.999 1.000 0.000 0.001 0.000
## 62 D80 0.999 1.000 0.000 0.001 0.000
## 63 D81 0.999 1.000 0.000 0.001 0.000
## 64 D82 0.998 1.000 0.000 0.001 0.000
## 65 D83 0.999 1.000 0.000 0.001 0.000
## 66 D84 0.998 1.000 0.000 0.001 0.000
## 67 D85 0.998 1.000 0.000 0.001 0.000
## 68 D86 0.998 1.000 0.000 0.001 0.000
## 69 D87 0.998 1.000 0.000 0.001 0.000
## 70 D88 0.998 1.000 0.000 0.001 0.000
## 71 D89 0.998 1.000 0.000 0.001 0.000
## 72 D90 0.997 1.000 0.000 0.001 0.000
## 73 D91 0.998 1.000 0.000 0.001 0.000
## 74 D92 0.997 1.000 0.000 0.001 0.000
## 75 D93 0.998 1.000 0.000 0.001 0.000
## 76 D94 0.998 1.000 0.000 0.001 0.000
## 77 J14 0.000 0.001 0.998 1.000 0.000
## 78 J15 0.000 0.001 0.997 1.000 0.000
## 79 J16 0.000 0.001 0.998 1.000 0.000
## 80 J17 0.000 0.001 0.999 1.000 0.000
## 81 J18 0.000 0.024 0.975 1.000 0.000
## 82 J19 0.000 0.001 0.999 1.000 0.000
## 83 J20 0.000 0.001 0.999 1.000 0.000
## 84 J21 0.000 0.001 0.999 1.000 0.000
## 85 J22 0.000 0.001 0.999 1.000 0.000
## 86 J23 0.000 0.001 0.998 1.000 0.000
## 87 J24 0.000 0.001 0.999 1.000 0.000
## 88 J25 0.000 0.001 0.999 1.000 0.000
## 89 J26 0.000 0.001 0.999 1.000 0.000
## 90 J27 0.000 0.001 0.999 1.000 0.000
## 91 K1 0.000 0.001 0.999 1.000 0.000
## 92 K2 0.000 0.001 0.999 1.000 0.000
## 93 K3 0.000 0.001 0.999 1.000 0.000
## 94 K32 0.000 0.001 0.999 1.000 0.000
## 95 K33 0.000 0.001 0.999 1.000 0.000
## 96 K34 0.000 0.001 0.998 1.000 0.000
## 97 K4 0.000 0.001 0.999 1.000 0.000
## 98 K5 0.000 0.001 0.998 1.000 0.000
## 99 L1 0.000 0.001 0.998 1.000 0.000
## 100 L10 0.000 0.001 0.991 1.000 0.000
## 101 L11 0.000 0.001 0.999 1.000 0.000
## 102 L12 0.000 0.001 0.998 1.000 0.000
## 103 L13 0.000 0.001 0.998 1.000 0.000
## 104 L14 0.000 0.001 0.999 1.000 0.000
## 105 L15 0.000 0.001 0.998 1.000 0.000
## 106 L16 0.000 0.001 0.999 1.000 0.000
## 107 L17 0.000 0.001 0.999 1.000 0.000
## 108 L18_rep 0.000 0.001 0.998 1.000 0.000
## 109 L19 0.000 0.001 0.999 1.000 0.000
## 110 L2 0.000 0.001 0.998 1.000 0.000
## 111 L20 0.000 0.001 0.999 1.000 0.000
## 112 L21 0.000 0.001 0.998 1.000 0.000
## 113 L22 0.000 0.001 0.999 1.000 0.000
## 114 L23 0.000 0.001 0.999 1.000 0.000
## 115 L24 0.000 0.001 0.998 1.000 0.000
## 116 L26 0.000 0.001 0.998 1.000 0.000
## 117 L27 0.000 0.001 0.998 1.000 0.000
## 118 L28 0.000 0.001 0.999 1.000 0.000
## 119 L3 0.000 0.001 0.998 1.000 0.000
## 120 L4 0.000 0.001 0.999 1.000 0.000
## 121 L5 0.000 0.001 0.998 1.000 0.000
## 122 L6 0.000 0.001 0.999 1.000 0.000
## 123 L7 0.000 0.001 0.995 1.000 0.000
## 124 L8 0.000 0.001 0.999 1.000 0.000
## 125 L9_rep 0.000 0.001 0.998 1.000 0.000
## 126 N1 0.000 0.001 0.995 1.000 0.000
## 127 N10 0.000 0.002 0.981 1.000 0.000
## 128 N11 0.000 0.001 0.997 1.000 0.000
## 129 N12 0.000 0.001 0.999 1.000 0.000
## 130 N13 0.000 0.001 0.999 1.000 0.000
## 131 N14 0.000 0.001 0.999 1.000 0.000
## 132 N15 0.000 0.001 0.999 1.000 0.000
## 133 N16 0.000 0.001 0.998 1.000 0.000
## 134 N17 0.000 0.001 0.999 1.000 0.000
## 135 N18 0.000 0.002 0.973 1.000 0.000
## 136 N19 0.000 0.001 0.998 1.000 0.000
## 137 N20 0.000 0.001 0.999 1.000 0.000
## 138 N21 0.000 0.001 0.998 1.000 0.000
## 139 N24 0.000 0.001 0.998 1.000 0.000
## 140 N25 0.000 0.001 0.998 1.000 0.000
## 141 N26 0.000 0.001 0.972 1.000 0.000
## 142 N27 0.000 0.001 0.987 1.000 0.000
## 143 N28 0.000 0.001 0.999 1.000 0.000
## 144 N3 0.000 0.002 0.987 1.000 0.000
## 145 N4 0.000 0.001 0.999 1.000 0.000
## 146 N5 0.000 0.001 0.998 1.000 0.000
## 147 N6 0.000 0.001 0.999 1.000 0.000
## 148 N8 0.000 0.001 0.998 1.000 0.000
## 149 N9 0.000 0.001 0.998 1.000 0.000
## 150 R10 0.000 0.001 0.999 1.000 0.000
## 151 R11 0.000 0.001 0.996 1.000 0.000
## 152 R12 0.000 0.001 0.998 1.000 0.000
## 153 R13 0.000 0.001 0.997 1.000 0.000
## 154 R14 0.000 0.001 0.999 1.000 0.000
## 155 R15 0.000 0.001 0.998 1.000 0.000
## 156 R16 0.000 0.001 0.999 1.000 0.000
## 157 R17 0.000 0.001 0.999 1.000 0.000
## 158 R18 0.000 0.001 0.999 1.000 0.000
## 159 R19 0.000 0.001 0.979 1.000 0.000
## 160 R2 0.000 0.001 0.998 1.000 0.000
## 161 R20 0.000 0.001 0.998 1.000 0.000
## 162 R3 0.000 0.001 0.998 1.000 0.000
## 163 R5 0.000 0.001 0.999 1.000 0.000
## 164 R6 0.000 0.001 0.999 1.000 0.000
## 165 R7 0.000 0.001 0.999 1.000 0.000
## 166 R9 0.000 0.001 0.999 1.000 0.000
## 167 S1 0.000 0.001 0.999 1.000 0.000
## 168 S11_rep 0.000 0.001 0.998 1.000 0.000
## 169 S12 0.000 0.001 0.998 1.000 0.000
## 170 S13 0.000 0.001 0.998 1.000 0.000
## 171 S15 0.000 0.001 0.999 1.000 0.000
## 172 S16 0.000 0.001 0.999 1.000 0.000
## 173 S17 0.000 0.001 0.999 1.000 0.000
## 174 S18 0.000 0.001 0.998 1.000 0.000
## 175 S19 0.000 0.001 0.999 1.000 0.000
## 176 S2 0.000 0.001 0.999 1.000 0.000
## 177 S20 0.000 0.001 0.998 1.000 0.000
## 178 S21 0.000 0.001 0.999 1.000 0.000
## 179 S22 0.000 0.001 0.998 1.000 0.000
## 180 S23 0.000 0.001 0.999 1.000 0.000
## 181 S24 0.000 0.001 0.999 1.000 0.000
## 182 S3 0.000 0.001 0.993 1.000 0.000
## 183 S4 0.000 0.001 0.999 1.000 0.000
## 184 S6 0.000 0.001 0.999 1.000 0.000
## 185 S7 0.000 0.001 0.998 1.000 0.000
## 186 S8 0.000 0.001 0.999 1.000 0.000
## 187 S9 0.000 0.001 0.998 1.000 0.000
## carinulata_H carinata_L carinata_H elongata sublineata carinulata carinata
## 1 0.112 0 0.001 0.000 0.912 0.088 0.000
## 2 0.089 0 0.001 0.931 0.000 0.068 0.000
## 3 0.046 0 0.001 0.970 0.000 0.030 0.000
## 4 0.001 0 0.001 1.000 0.000 0.000 0.000
## 5 0.001 0 0.001 1.000 0.000 0.000 0.000
## 6 0.001 0 0.001 1.000 0.000 0.000 0.000
## 7 0.001 0 0.001 1.000 0.000 0.000 0.000
## 8 0.001 0 0.005 0.998 0.001 0.000 0.001
## 9 0.001 0 0.001 1.000 0.000 0.000 0.000
## 10 0.001 0 0.001 1.000 0.000 0.000 0.000
## 11 0.001 0 0.001 1.000 0.000 0.000 0.000
## 12 0.001 0 0.001 1.000 0.000 0.000 0.000
## 13 0.001 0 0.001 1.000 0.000 0.000 0.000
## 14 0.001 0 0.001 1.000 0.000 0.000 0.000
## 15 0.001 0 0.001 1.000 0.000 0.000 0.000
## 16 0.001 0 0.001 1.000 0.000 0.000 0.000
## 17 0.001 0 0.001 1.000 0.000 0.000 0.000
## 18 0.001 0 0.001 1.000 0.000 0.000 0.000
## 19 0.001 0 0.001 1.000 0.000 0.000 0.000
## 20 0.001 0 0.001 1.000 0.000 0.000 0.000
## 21 0.001 0 0.001 1.000 0.000 0.000 0.000
## 22 0.001 0 0.001 1.000 0.000 0.000 0.000
## 23 0.001 0 0.001 1.000 0.000 0.000 0.000
## 24 0.001 0 0.001 1.000 0.000 0.000 0.000
## 25 0.001 0 0.001 1.000 0.000 0.000 0.000
## 26 0.001 0 0.001 1.000 0.000 0.000 0.000
## 27 0.001 0 0.001 1.000 0.000 0.000 0.000
## 28 0.001 0 0.001 1.000 0.000 0.000 0.000
## 29 0.001 0 0.001 1.000 0.000 0.000 0.000
## 30 0.001 0 0.001 1.000 0.000 0.000 0.000
## 31 0.001 0 0.001 1.000 0.000 0.000 0.000
## 32 0.001 0 0.001 1.000 0.000 0.000 0.000
## 33 0.001 0 0.001 1.000 0.000 0.000 0.000
## 34 0.001 0 0.001 1.000 0.000 0.000 0.000
## 35 0.001 0 0.001 1.000 0.000 0.000 0.000
## 36 0.001 0 0.001 1.000 0.000 0.000 0.000
## 37 0.001 0 0.001 1.000 0.000 0.000 0.000
## 38 0.001 0 0.001 1.000 0.000 0.000 0.000
## 39 0.001 0 0.001 1.000 0.000 0.000 0.000
## 40 0.001 0 0.001 1.000 0.000 0.000 0.000
## 41 0.001 0 0.001 1.000 0.000 0.000 0.000
## 42 0.001 0 0.001 1.000 0.000 0.000 0.000
## 43 0.001 0 0.001 1.000 0.000 0.000 0.000
## 44 0.001 0 0.001 1.000 0.000 0.000 0.000
## 45 0.001 0 0.001 1.000 0.000 0.000 0.000
## 46 0.001 0 0.001 1.000 0.000 0.000 0.000
## 47 0.001 0 0.001 1.000 0.000 0.000 0.000
## 48 0.001 0 0.001 1.000 0.000 0.000 0.000
## 49 0.001 0 0.001 1.000 0.000 0.000 0.000
## 50 0.001 0 0.001 1.000 0.000 0.000 0.000
## 51 0.029 0 0.001 0.989 0.002 0.010 0.000
## 52 0.001 0 0.001 1.000 0.000 0.000 0.000
## 53 0.012 0 0.001 0.998 0.000 0.002 0.000
## 54 0.001 0 0.001 1.000 0.000 0.000 0.000
## 55 0.001 0 0.001 1.000 0.000 0.000 0.000
## 56 0.001 0 0.001 1.000 0.000 0.000 0.000
## 57 0.001 0 0.001 1.000 0.000 0.000 0.000
## 58 0.001 0 0.001 1.000 0.000 0.000 0.000
## 59 0.001 0 0.001 1.000 0.000 0.000 0.000
## 60 0.001 0 0.001 1.000 0.000 0.000 0.000
## 61 0.001 0 0.001 1.000 0.000 0.000 0.000
## 62 0.001 0 0.001 1.000 0.000 0.000 0.000
## 63 0.001 0 0.001 1.000 0.000 0.000 0.000
## 64 0.001 0 0.001 1.000 0.000 0.000 0.000
## 65 0.001 0 0.001 1.000 0.000 0.000 0.000
## 66 0.001 0 0.001 1.000 0.000 0.000 0.000
## 67 0.001 0 0.001 1.000 0.000 0.000 0.000
## 68 0.001 0 0.001 1.000 0.000 0.000 0.000
## 69 0.001 0 0.001 1.000 0.000 0.000 0.000
## 70 0.001 0 0.001 1.000 0.000 0.000 0.000
## 71 0.001 0 0.001 1.000 0.000 0.000 0.000
## 72 0.001 0 0.001 1.000 0.000 0.000 0.000
## 73 0.001 0 0.001 1.000 0.000 0.000 0.000
## 74 0.001 0 0.001 1.000 0.000 0.000 0.000
## 75 0.001 0 0.001 1.000 0.000 0.000 0.000
## 76 0.001 0 0.001 1.000 0.000 0.000 0.000
## 77 0.001 0 0.001 0.000 1.000 0.000 0.000
## 78 0.001 0 0.002 0.000 1.000 0.000 0.000
## 79 0.001 0 0.001 0.000 1.000 0.000 0.000
## 80 0.001 0 0.001 0.000 1.000 0.000 0.000
## 81 0.001 0 0.001 0.012 0.988 0.000 0.000
## 82 0.001 0 0.001 0.000 1.000 0.000 0.000
## 83 0.001 0 0.001 0.000 1.000 0.000 0.000
## 84 0.001 0 0.001 0.000 1.000 0.000 0.000
## 85 0.001 0 0.001 0.000 1.000 0.000 0.000
## 86 0.001 0 0.001 0.000 1.000 0.000 0.000
## 87 0.001 0 0.001 0.000 1.000 0.000 0.000
## 88 0.001 0 0.001 0.000 1.000 0.000 0.000
## 89 0.001 0 0.001 0.000 1.000 0.000 0.000
## 90 0.001 0 0.001 0.000 1.000 0.000 0.000
## 91 0.001 0 0.001 0.000 1.000 0.000 0.000
## 92 0.001 0 0.001 0.000 1.000 0.000 0.000
## 93 0.001 0 0.001 0.000 1.000 0.000 0.000
## 94 0.001 0 0.001 0.000 1.000 0.000 0.000
## 95 0.001 0 0.001 0.000 1.000 0.000 0.000
## 96 0.001 0 0.001 0.000 1.000 0.000 0.000
## 97 0.001 0 0.001 0.000 1.000 0.000 0.000
## 98 0.001 0 0.001 0.000 1.000 0.000 0.000
## 99 0.001 0 0.001 0.000 1.000 0.000 0.000
## 100 0.004 0 0.005 0.000 0.999 0.000 0.001
## 101 0.001 0 0.001 0.000 1.000 0.000 0.000
## 102 0.001 0 0.001 0.000 1.000 0.000 0.000
## 103 0.001 0 0.001 0.000 1.000 0.000 0.000
## 104 0.001 0 0.001 0.000 1.000 0.000 0.000
## 105 0.001 0 0.001 0.000 1.000 0.000 0.000
## 106 0.001 0 0.001 0.000 1.000 0.000 0.000
## 107 0.001 0 0.001 0.000 1.000 0.000 0.000
## 108 0.001 0 0.001 0.000 1.000 0.000 0.000
## 109 0.001 0 0.001 0.000 1.000 0.000 0.000
## 110 0.001 0 0.001 0.000 1.000 0.000 0.000
## 111 0.001 0 0.001 0.000 1.000 0.000 0.000
## 112 0.001 0 0.001 0.000 1.000 0.000 0.000
## 113 0.001 0 0.001 0.000 1.000 0.000 0.000
## 114 0.001 0 0.001 0.000 1.000 0.000 0.000
## 115 0.001 0 0.001 0.000 1.000 0.000 0.000
## 116 0.001 0 0.001 0.000 1.000 0.000 0.000
## 117 0.001 0 0.001 0.000 1.000 0.000 0.000
## 118 0.001 0 0.001 0.000 1.000 0.000 0.000
## 119 0.001 0 0.001 0.000 1.000 0.000 0.000
## 120 0.001 0 0.001 0.000 1.000 0.000 0.000
## 121 0.001 0 0.001 0.000 1.000 0.000 0.000
## 122 0.001 0 0.001 0.000 1.000 0.000 0.000
## 123 0.001 0 0.005 0.000 0.999 0.000 0.001
## 124 0.001 0 0.001 0.000 1.000 0.000 0.000
## 125 0.001 0 0.001 0.000 1.000 0.000 0.000
## 126 0.001 0 0.004 0.000 0.999 0.000 0.001
## 127 0.018 0 0.002 0.000 0.994 0.006 0.000
## 128 0.002 0 0.001 0.000 1.000 0.000 0.000
## 129 0.001 0 0.001 0.000 1.000 0.000 0.000
## 130 0.001 0 0.001 0.000 1.000 0.000 0.000
## 131 0.001 0 0.001 0.000 1.000 0.000 0.000
## 132 0.001 0 0.001 0.000 1.000 0.000 0.000
## 133 0.001 0 0.001 0.000 1.000 0.000 0.000
## 134 0.001 0 0.001 0.000 1.000 0.000 0.000
## 135 0.023 0 0.017 0.000 0.989 0.008 0.002
## 136 0.001 0 0.001 0.000 1.000 0.000 0.000
## 137 0.001 0 0.001 0.000 1.000 0.000 0.000
## 138 0.001 0 0.001 0.000 1.000 0.000 0.000
## 139 0.001 0 0.001 0.000 1.000 0.000 0.000
## 140 0.001 0 0.001 0.000 1.000 0.000 0.000
## 141 0.027 0 0.004 0.000 0.986 0.013 0.001
## 142 0.013 0 0.001 0.000 0.997 0.002 0.000
## 143 0.001 0 0.001 0.000 1.000 0.000 0.000
## 144 0.012 0 0.001 0.000 0.998 0.002 0.000
## 145 0.001 0 0.001 0.000 1.000 0.000 0.000
## 146 0.001 0 0.001 0.000 1.000 0.000 0.000
## 147 0.001 0 0.001 0.000 1.000 0.000 0.000
## 148 0.001 0 0.001 0.000 1.000 0.000 0.000
## 149 0.001 0 0.001 0.000 1.000 0.000 0.000
## 150 0.001 0 0.001 0.000 1.000 0.000 0.000
## 151 0.001 0 0.003 0.000 0.999 0.000 0.000
## 152 0.001 0 0.001 0.000 1.000 0.000 0.000
## 153 0.001 0 0.002 0.000 0.999 0.000 0.000
## 154 0.001 0 0.001 0.000 1.000 0.000 0.000
## 155 0.001 0 0.001 0.000 1.000 0.000 0.000
## 156 0.001 0 0.001 0.000 1.000 0.000 0.000
## 157 0.001 0 0.001 0.000 1.000 0.000 0.000
## 158 0.001 0 0.001 0.000 1.000 0.000 0.000
## 159 0.016 0 0.013 0.000 0.994 0.004 0.002
## 160 0.001 0 0.001 0.000 1.000 0.000 0.000
## 161 0.001 0 0.001 0.000 1.000 0.000 0.000
## 162 0.001 0 0.001 0.000 1.000 0.000 0.000
## 163 0.001 0 0.001 0.000 1.000 0.000 0.000
## 164 0.001 0 0.001 0.000 1.000 0.000 0.000
## 165 0.001 0 0.001 0.000 1.000 0.000 0.000
## 166 0.001 0 0.001 0.000 1.000 0.000 0.000
## 167 0.001 0 0.001 0.000 1.000 0.000 0.000
## 168 0.001 0 0.001 0.000 1.000 0.000 0.000
## 169 0.001 0 0.001 0.000 1.000 0.000 0.000
## 170 0.001 0 0.001 0.000 1.000 0.000 0.000
## 171 0.001 0 0.001 0.000 1.000 0.000 0.000
## 172 0.001 0 0.001 0.000 1.000 0.000 0.000
## 173 0.001 0 0.001 0.000 1.000 0.000 0.000
## 174 0.001 0 0.001 0.000 1.000 0.000 0.000
## 175 0.001 0 0.001 0.000 1.000 0.000 0.000
## 176 0.001 0 0.001 0.000 1.000 0.000 0.000
## 177 0.001 0 0.001 0.000 1.000 0.000 0.000
## 178 0.001 0 0.001 0.000 1.000 0.000 0.000
## 179 0.001 0 0.001 0.000 1.000 0.000 0.000
## 180 0.001 0 0.001 0.000 1.000 0.000 0.000
## 181 0.001 0 0.001 0.000 1.000 0.000 0.000
## 182 0.002 0 0.003 0.000 0.999 0.000 0.001
## 183 0.001 0 0.001 0.000 1.000 0.000 0.000
## 184 0.001 0 0.001 0.000 1.000 0.000 0.000
## 185 0.001 0 0.001 0.000 1.000 0.000 0.000
## 186 0.001 0 0.001 0.000 1.000 0.000 0.000
## 187 0.001 0 0.001 0.000 1.000 0.000 0.000
## popID
## 1 47TX
## 2 44GR
## 3 41CR
## 4 38CR
## 5 41CR
## 6 41CR
## 7 41CR
## 8 43GR
## 9 43GR
## 10 43GR
## 11 38CR
## 12 43GR
## 13 44GR
## 14 44GR
## 15 44GR
## 16 44GR
## 17 38CR
## 18 43GR
## 19 43GR
## 20 43GR
## 21 44GR
## 22 44GR
## 23 44GR
## 24 38CR
## 25 44GR
## 26 41CR
## 27 41CR
## 28 41CR
## 29 41CR
## 30 37CR
## 31 37CR
## 32 37CR
## 33 38CR
## 34 39CR
## 35 38CR
## 36 38CR
## 37 38CR
## 38 39CR
## 39 39CR
## 40 39CR
## 41 39CR
## 42 41CR
## 43 41CR
## 44 41CR
## 45 39CR
## 46 41CR
## 47 43GR
## 48 43GR
## 49 43GR
## 50 43GR
## 51 44GR
## 52 39CR
## 53 37CR
## 54 37CR
## 55 38CR
## 56 38CR
## 57 38CR
## 58 39CR
## 59 39CR
## 60 39CR
## 61 39CR
## 62 39CR
## 63 41CR
## 64 41CR
## 65 41CR
## 66 41CR
## 67 43GR
## 68 43GR
## 69 43GR
## 70 43GR
## 71 44GR
## 72 44GR
## 73 44GR
## 74 44GR
## 75 44GR
## 76 44GR
## 77 51TX
## 78 51TX
## 79 51TX
## 80 51TX
## 81 51TX
## 82 51TX
## 83 51TX
## 84 51TX
## 85 51TX
## 86 51TX
## 87 51TX
## 88 51TX
## 89 51TX
## 90 51TX
## 91 50TX
## 92 50TX
## 93 50TX
## 94 50TX
## 95 50TX
## 96 50TX
## 97 50TX
## 98 50TX
## 99 49TX
## 100 49TX
## 101 49TX
## 102 49TX
## 103 49TX
## 104 49TX
## 105 49TX
## 106 49TX
## 107 49TX
## 108 49TX
## 109 49TX
## 110 49TX
## 111 49TX
## 112 49TX
## 113 49TX
## 114 49TX
## 115 49TX
## 116 49TX
## 117 49TX
## 118 49TX
## 119 49TX
## 120 49TX
## 121 49TX
## 122 49TX
## 123 49TX
## 124 49TX
## 125 49TX
## 126 47TX
## 127 47TX
## 128 47TX
## 129 47TX
## 130 47TX
## 131 47TX
## 132 47TX
## 133 47TX
## 134 47TX
## 135 47TX
## 136 47TX
## 137 47TX
## 138 47TX
## 139 47TX
## 140 47TX
## 141 47TX
## 142 47TX
## 143 47TX
## 144 47TX
## 145 47TX
## 146 47TX
## 147 47TX
## 148 47TX
## 149 47TX
## 150 SUB_LAB
## 151 SUB_LAB
## 152 SUB_LAB
## 153 SUB_LAB
## 154 SUB_LAB
## 155 SUB_LAB
## 156 SUB_LAB
## 157 SUB_LAB
## 158 SUB_LAB
## 159 SUB_LAB
## 160 SUB_LAB
## 161 SUB_LAB
## 162 SUB_LAB
## 163 SUB_LAB
## 164 SUB_LAB
## 165 SUB_LAB
## 166 SUB_LAB
## 167 52TX
## 168 52TX
## 169 52TX
## 170 52TX
## 171 52TX
## 172 52TX
## 173 52TX
## 174 52TX
## 175 52TX
## 176 52TX
## 177 52TX
## 178 52TX
## 179 52TX
## 180 52TX
## 181 52TX
## 182 52TX
## 183 52TX
## 184 52TX
## 185 52TX
## 186 52TX
## 187 52TXfind_threshold <- ci %>%
semi_join(diagnostic) %>%
arrange(desc(carinulata_L))## Joining, by = "popID"t <- find_threshold[1,"carinulata"]
t## [1] 0.088t_low <- find_threshold[1,"carinulata_L"]
t_low## [1] 0.066Each hybrid type
How many hybrids with carinulata?
carinu_hyb <- ci %>%
dplyr::filter(carinulata_L > t_low ) %>%
dplyr::filter(elongata_L > t_low |
sublineata_L > t_low |
carinata_L > t_low) %>%
arrange(desc(carinulata_L))
carinu_hyb## sampID elongata_L elongata_H sublineata_L sublineata_H carinulata_L
## 1 G5_rep 0 0.001 0.031 0.116 0.633
## 2 G48_rep 0 0.001 0.890 0.932 0.068
## carinulata_H carinata_L carinata_H elongata sublineata carinulata carinata
## 1 0.695 0.216 0.312 0 0.071 0.664 0.265
## 2 0.109 0.000 0.001 0 0.912 0.088 0.000
## popID
## 1 28NM
## 2 6KSHow many tri-specific hybrids?
ci %>%
dplyr::filter(carinata_L > t_low & elongata_L > t_low & sublineata_L > t_low) ## [1] sampID elongata_L elongata_H sublineata_L sublineata_H
## [6] carinulata_L carinulata_H carinata_L carinata_H elongata
## [11] sublineata carinulata carinata popID
## <0 rows> (or 0-length row.names)## If we count less-conservatively:
ci %>%
dplyr::filter(carinata_L >.01 & elongata_L> .01 & sublineata_L > .01) %>%
arrange(desc(elongata_L))## sampID elongata_L elongata_H sublineata_L sublineata_H carinulata_L
## 1 H10 0.024 0.079 0.140 0.258 0.015
## 2 F84 0.018 0.052 0.867 0.940 0.000
## 3 G2 0.018 0.051 0.111 0.202 0.000
## carinulata_H carinata_L carinata_H elongata sublineata carinulata carinata
## 1 0.056 0.655 0.777 0.050 0.198 0.035 0.717
## 2 0.001 0.024 0.100 0.034 0.905 0.000 0.060
## 3 0.001 0.762 0.856 0.034 0.155 0.000 0.811
## popID
## 1 28NM
## 2 26NM
## 3 15TXHow many hybrids with sublineata?
ci %>%
dplyr::filter(sublineata_L > t_low) %>%
dplyr::filter(carinulata_L > t_low |
elongata_L > t_low | carinata_L > t_low) %>%
arrange(desc(sublineata_L)) ## sampID elongata_L elongata_H sublineata_L sublineata_H carinulata_L
## 1 G48_rep 0.000 0.001 0.890 0.932 0.068
## 2 G10 0.088 0.133 0.866 0.912 0.000
## 3 G50 0.092 0.140 0.860 0.908 0.000
## 4 F95 0.000 0.001 0.844 0.928 0.000
## 5 G11 0.000 0.001 0.818 0.912 0.000
## 6 G6 0.000 0.001 0.815 0.897 0.000
## 7 H21 0.000 0.001 0.797 0.887 0.000
## 8 H71 0.000 0.036 0.796 0.885 0.000
## 9 H23 0.000 0.037 0.793 0.878 0.000
## 10 G63 0.000 0.019 0.762 0.869 0.000
## 11 H26 0.198 0.263 0.736 0.802 0.000
## 12 F94 0.002 0.051 0.727 0.838 0.000
## 13 H12 0.000 0.001 0.727 0.828 0.000
## 14 G9 0.000 0.001 0.721 0.820 0.000
## 15 G52 0.000 0.021 0.667 0.771 0.000
## 16 G83 0.000 0.002 0.642 0.802 0.000
## 17 F96 0.000 0.003 0.619 0.739 0.000
## 18 G19 0.000 0.001 0.572 0.704 0.000
## 19 H67 0.000 0.001 0.472 0.599 0.000
## 20 H68 0.000 0.001 0.397 0.521 0.000
## 21 F92 0.000 0.001 0.380 0.511 0.000
## 22 H22 0.000 0.007 0.378 0.503 0.000
## 23 H9 0.000 0.002 0.348 0.466 0.000
## 24 H64 0.000 0.001 0.291 0.412 0.000
## 25 G38 0.000 0.001 0.271 0.393 0.000
## 26 G78 0.000 0.001 0.154 0.262 0.000
## 27 H10 0.024 0.079 0.140 0.258 0.015
## 28 G2 0.018 0.051 0.111 0.202 0.000
## carinulata_H carinata_L carinata_H elongata sublineata carinulata carinata
## 1 0.109 0.000 0.001 0.000 0.912 0.088 0.000
## 2 0.001 0.000 0.003 0.110 0.890 0.000 0.001
## 3 0.001 0.000 0.001 0.115 0.885 0.000 0.000
## 4 0.001 0.071 0.156 0.000 0.888 0.000 0.111
## 5 0.001 0.088 0.182 0.000 0.867 0.000 0.133
## 6 0.001 0.103 0.185 0.000 0.858 0.000 0.142
## 7 0.001 0.113 0.202 0.000 0.844 0.000 0.156
## 8 0.001 0.096 0.190 0.017 0.843 0.000 0.140
## 9 0.001 0.100 0.189 0.020 0.837 0.000 0.143
## 10 0.004 0.125 0.235 0.003 0.817 0.000 0.179
## 11 0.001 0.000 0.001 0.230 0.770 0.000 0.000
## 12 0.001 0.130 0.247 0.028 0.785 0.000 0.187
## 13 0.001 0.172 0.272 0.000 0.779 0.000 0.221
## 14 0.001 0.180 0.279 0.000 0.772 0.000 0.228
## 15 0.001 0.223 0.330 0.004 0.720 0.000 0.276
## 16 0.001 0.198 0.357 0.000 0.725 0.000 0.274
## 17 0.001 0.260 0.381 0.001 0.680 0.000 0.320
## 18 0.001 0.296 0.428 0.000 0.639 0.000 0.360
## 19 0.001 0.401 0.527 0.000 0.536 0.000 0.464
## 20 0.001 0.479 0.603 0.000 0.459 0.000 0.541
## 21 0.001 0.489 0.620 0.000 0.445 0.000 0.555
## 22 0.001 0.496 0.621 0.001 0.440 0.000 0.559
## 23 0.001 0.533 0.651 0.000 0.407 0.000 0.593
## 24 0.001 0.588 0.709 0.000 0.350 0.000 0.650
## 25 0.001 0.607 0.729 0.000 0.331 0.000 0.668
## 26 0.001 0.737 0.846 0.000 0.207 0.000 0.793
## 27 0.056 0.655 0.777 0.050 0.198 0.035 0.717
## 28 0.001 0.762 0.856 0.034 0.155 0.000 0.811
## popID
## 1 6KS
## 2 26NM
## 3 20NM
## 4 26NM
## 5 26NM
## 6 28NM
## 7 26NM
## 8 26NM
## 9 26NM
## 10 20NM
## 11 21NM
## 12 26NM
## 13 26NM
## 14 26NM
## 15 20NM
## 16 21NM
## 17 26NM
## 18 27NM
## 19 28NM
## 20 28NM
## 21 28NM
## 22 26NM
## 23 28NM
## 24 18TX
## 25 20NM
## 26 16TX
## 27 28NM
## 28 15TXci %>%
dplyr::filter(sublineata_L > t) %>%
dplyr::filter(carinulata_L > t_low |
elongata_L > t_low | carinata_L > t_low) %>%
arrange(desc(sublineata_L)) %>%
count()## n
## 1 28D. elongata hybrids
ci %>%
dplyr::filter(elongata_L > t_low) %>%
dplyr::filter(carinulata_L > t_low |
sublineata_L > t_low | carinata_L > t_low) %>%
arrange(desc(elongata_L))## sampID elongata_L elongata_H sublineata_L sublineata_H carinulata_L
## 1 E29 0.221 0.305 0.000 0.002 0
## 2 H26 0.198 0.263 0.736 0.802 0
## 3 G50 0.092 0.140 0.860 0.908 0
## 4 G10 0.088 0.133 0.866 0.912 0
## 5 G87 0.082 0.137 0.000 0.005 0
## 6 G01_rep 0.078 0.129 0.000 0.001 0
## carinulata_H carinata_L carinata_H elongata sublineata carinulata carinata
## 1 0.001 0.694 0.779 0.262 0.000 0 0.737
## 2 0.001 0.000 0.001 0.230 0.770 0 0.000
## 3 0.001 0.000 0.001 0.115 0.885 0 0.000
## 4 0.001 0.000 0.003 0.110 0.890 0 0.001
## 5 0.001 0.862 0.917 0.109 0.001 0 0.890
## 6 0.001 0.870 0.922 0.103 0.000 0 0.897
## popID
## 1 15TX
## 2 21NM
## 3 20NM
## 4 26NM
## 5 15TX
## 6 15TXci %>%
dplyr::filter(elongata_L > t_low) %>%
dplyr::filter(carinulata_L > t_low |
sublineata_L > t_low | carinata_L > t_low) %>%
arrange(desc(elongata_L)) %>%
count()## n
## 1 6Specific bi-parental types
How many carinata-sublineata hybrids? How many carinata-elongata hybrids? How many elong-sublineata hybrids?
car_sub_hyb <- ci %>%
dplyr::filter(carinata_L >t_low & sublineata_L > t_low & elongata_L < t_low)
car_sub_hyb## sampID elongata_L elongata_H sublineata_L sublineata_H carinulata_L
## 1 F92 0.000 0.001 0.380 0.511 0.000
## 2 F94 0.002 0.051 0.727 0.838 0.000
## 3 F95 0.000 0.001 0.844 0.928 0.000
## 4 F96 0.000 0.003 0.619 0.739 0.000
## 5 G11 0.000 0.001 0.818 0.912 0.000
## 6 G19 0.000 0.001 0.572 0.704 0.000
## 7 G2 0.018 0.051 0.111 0.202 0.000
## 8 G38 0.000 0.001 0.271 0.393 0.000
## 9 G52 0.000 0.021 0.667 0.771 0.000
## 10 G6 0.000 0.001 0.815 0.897 0.000
## 11 G63 0.000 0.019 0.762 0.869 0.000
## 12 G78 0.000 0.001 0.154 0.262 0.000
## 13 G83 0.000 0.002 0.642 0.802 0.000
## 14 G9 0.000 0.001 0.721 0.820 0.000
## 15 H10 0.024 0.079 0.140 0.258 0.015
## 16 H12 0.000 0.001 0.727 0.828 0.000
## 17 H21 0.000 0.001 0.797 0.887 0.000
## 18 H22 0.000 0.007 0.378 0.503 0.000
## 19 H23 0.000 0.037 0.793 0.878 0.000
## 20 H64 0.000 0.001 0.291 0.412 0.000
## 21 H67 0.000 0.001 0.472 0.599 0.000
## 22 H68 0.000 0.001 0.397 0.521 0.000
## 23 H71 0.000 0.036 0.796 0.885 0.000
## 24 H9 0.000 0.002 0.348 0.466 0.000
## carinulata_H carinata_L carinata_H elongata sublineata carinulata carinata
## 1 0.001 0.489 0.620 0.000 0.445 0.000 0.555
## 2 0.001 0.130 0.247 0.028 0.785 0.000 0.187
## 3 0.001 0.071 0.156 0.000 0.888 0.000 0.111
## 4 0.001 0.260 0.381 0.001 0.680 0.000 0.320
## 5 0.001 0.088 0.182 0.000 0.867 0.000 0.133
## 6 0.001 0.296 0.428 0.000 0.639 0.000 0.360
## 7 0.001 0.762 0.856 0.034 0.155 0.000 0.811
## 8 0.001 0.607 0.729 0.000 0.331 0.000 0.668
## 9 0.001 0.223 0.330 0.004 0.720 0.000 0.276
## 10 0.001 0.103 0.185 0.000 0.858 0.000 0.142
## 11 0.004 0.125 0.235 0.003 0.817 0.000 0.179
## 12 0.001 0.737 0.846 0.000 0.207 0.000 0.793
## 13 0.001 0.198 0.357 0.000 0.725 0.000 0.274
## 14 0.001 0.180 0.279 0.000 0.772 0.000 0.228
## 15 0.056 0.655 0.777 0.050 0.198 0.035 0.717
## 16 0.001 0.172 0.272 0.000 0.779 0.000 0.221
## 17 0.001 0.113 0.202 0.000 0.844 0.000 0.156
## 18 0.001 0.496 0.621 0.001 0.440 0.000 0.559
## 19 0.001 0.100 0.189 0.020 0.837 0.000 0.143
## 20 0.001 0.588 0.709 0.000 0.350 0.000 0.650
## 21 0.001 0.401 0.527 0.000 0.536 0.000 0.464
## 22 0.001 0.479 0.603 0.000 0.459 0.000 0.541
## 23 0.001 0.096 0.190 0.017 0.843 0.000 0.140
## 24 0.001 0.533 0.651 0.000 0.407 0.000 0.593
## popID
## 1 28NM
## 2 26NM
## 3 26NM
## 4 26NM
## 5 26NM
## 6 27NM
## 7 15TX
## 8 20NM
## 9 20NM
## 10 28NM
## 11 20NM
## 12 16TX
## 13 21NM
## 14 26NM
## 15 28NM
## 16 26NM
## 17 26NM
## 18 26NM
## 19 26NM
## 20 18TX
## 21 28NM
## 22 28NM
## 23 26NM
## 24 28NMcar_sub_hyb %>%
summarise(mean(sublineata), nrow(.))## mean(sublineata) nrow(.)
## 1 0.6075833 24car_elong_hyb <- ci %>%
dplyr::filter(carinata_L > t_low & elongata_L> t_low & sublineata_L < t_low) %>%
arrange(desc(elongata_L))
car_elong_hyb## sampID elongata_L elongata_H sublineata_L sublineata_H carinulata_L
## 1 E29 0.221 0.305 0 0.002 0
## 2 G87 0.082 0.137 0 0.005 0
## 3 G01_rep 0.078 0.129 0 0.001 0
## carinulata_H carinata_L carinata_H elongata sublineata carinulata carinata
## 1 0.001 0.694 0.779 0.262 0.000 0 0.737
## 2 0.001 0.862 0.917 0.109 0.001 0 0.890
## 3 0.001 0.870 0.922 0.103 0.000 0 0.897
## popID
## 1 15TX
## 2 15TX
## 3 15TXcar_elong_hyb %>%
summarise(mean(carinata), nrow(.))## mean(carinata) nrow(.)
## 1 0.8413333 3elong_sub_hyb <- ci %>%
dplyr::filter(carinata_L < t_low & sublineata_L > t_low & elongata_L > t_low)
elong_sub_hyb## sampID elongata_L elongata_H sublineata_L sublineata_H carinulata_L
## 1 G10 0.088 0.133 0.866 0.912 0
## 2 G50 0.092 0.140 0.860 0.908 0
## 3 H26 0.198 0.263 0.736 0.802 0
## carinulata_H carinata_L carinata_H elongata sublineata carinulata carinata
## 1 0.001 0 0.003 0.110 0.890 0 0.001
## 2 0.001 0 0.001 0.115 0.885 0 0.000
## 3 0.001 0 0.001 0.230 0.770 0 0.000
## popID
## 1 26NM
## 2 20NM
## 3 21NMelong_sub_hyb %>%
summarise(mean(sublineata), nrow(.))## mean(sublineata) nrow(.)
## 1 0.8483333 3Vizualize bi-parental crosses!
This could be better represented as a simplex (triangle)?
pCS <- ggplot(car_sub_hyb) +
geom_histogram(aes(carinata),
alpha = .5, fill = legend[which(legend$spp=="carinata"), "cbbPalette"], stat = "bin", binwidth = .05) +
geom_histogram(aes(sublineata),
alpha = .5, fill = legend[which(legend$spp=="sublineata"), "cbbPalette"], stat = "bin", binwidth = .05) +
ylim(c(0,8)) +
scale_x_continuous(name= expression( ~italic("D. carinata")~ "X" ~italic("D. sublineata")~ "q-values")) +
theme_bw()
pCS
pCE <-ggplot(car_elong_hyb) +
geom_histogram(aes(carinata),
alpha = .5, fill = legend[which(legend$spp=="carinata"), "cbbPalette"], stat = "bin", binwidth = .05) +
geom_histogram(aes(elongata),
alpha = .5, fill = legend[which(legend$spp=="elongata"), "cbbPalette"], stat = "bin", binwidth = .05) +
ylim(c(0,8)) +
scale_x_continuous(name= expression( ~italic("D. elongata")~ "X" ~italic(" D. carinata")~ "q-values")) +
theme_bw()
pCE
pES <-ggplot(elong_sub_hyb) +
geom_histogram(aes(elongata),
alpha = .5,fill = legend[which(legend$spp=="elongata"), "cbbPalette"], stat = "bin", binwidth = .05) +
geom_histogram(aes(sublineata),
alpha = .5, fill = legend[which(legend$spp=="sublineata"), "cbbPalette"], stat = "bin", binwidth = .05) +
ylim(c(0,8)) +
scale_x_continuous(name= expression( ~italic("D. elongata")~ "X" ~italic(" D. sublineata")~ "q-values")) +
theme_bw()
cowplot::plot_grid(pCS, pCE, pES, nrow = 1, labels = c('A','B','C'))
ggsave("hybrid_distributions.jpg", dpi = 300, height = 4, width = 10, units = "in")Which sites have individuals with admxiture?
sample_size <- samplesize %>%
select('popID', 'n')
hybrids <- rbind(elong_sub_hyb, car_elong_hyb, carinu_hyb, car_sub_hyb)
nrow(hybrids)## [1] 32no_hybrids <- anti_join(ci, hybrids)## Joining, by = c("sampID", "elongata_L", "elongata_H", "sublineata_L", "sublineata_H", "carinulata_L", "carinulata_H", "carinata_L", "carinata_H", "elongata", "sublineata", "carinulata", "carinata", "popID")pure_carinu <- dplyr::filter(no_hybrids, carinulata > .9)
nrow(pure_carinu)## [1] 179pure_sub <- dplyr::filter(no_hybrids, sublineata > .9)
nrow(pure_sub)## [1] 176pure_carinata <- dplyr::filter(no_hybrids, carinata > .9 )
nrow(pure_carinata)## [1] 69pure_elongata <- dplyr::filter(no_hybrids, elongata > .9 )
nrow(pure_elongata) ## [1] 93pure_parentals <- rbind(pure_carinu, pure_sub, pure_carinata, pure_elongata)
what_else <- anti_join( ci, pure_parentals) %>%
anti_join(., hybrids)## Joining, by = c("sampID", "elongata_L", "elongata_H", "sublineata_L", "sublineata_H", "carinulata_L", "carinulata_H", "carinata_L", "carinata_H", "elongata", "sublineata", "carinulata", "carinata", "popID")
## Joining, by = c("sampID", "elongata_L", "elongata_H", "sublineata_L", "sublineata_H", "carinulata_L", "carinulata_H", "carinata_L", "carinata_H", "elongata", "sublineata", "carinulata", "carinata", "popID") num_admixed_inds <- hybrids %>%
group_by(popID) %>%
count(popID)
colnames(num_admixed_inds)[2] <- 'num_admixed'
num_admixed_inds## # A tibble: 9 x 2
## # Groups: popID [9]
## popID num_admixed
## <chr> <int>
## 1 15TX 4
## 2 16TX 1
## 3 18TX 1
## 4 20NM 4
## 5 21NM 2
## 6 26NM 11
## 7 27NM 1
## 8 28NM 7
## 9 6KS 1merge(sample_size, num_admixed_inds) %>%
mutate(perc_admixed = num_admixed/n) %>%
arrange(desc(perc_admixed))## popID n num_admixed perc_admixed
## 1 26NM 18 11 0.61111111
## 2 28NM 17 7 0.41176471
## 3 15TX 16 4 0.25000000
## 4 20NM 17 4 0.23529412
## 5 21NM 14 2 0.14285714
## 6 6KS 10 1 0.10000000
## 7 18TX 17 1 0.05882353
## 8 16TX 18 1 0.05555556
## 9 27NM 20 1 0.05000000hybrids %>%
filter_at(vars(elongata:carinata),
~ . < .95) %>%
group_by(popID) %>%
summarize_at(vars(carinata:sublineata), funs(mean, var))## Warning: `funs()` was deprecated in dplyr 0.8.0.
## Please use a list of either functions or lambdas:
##
## # Simple named list:
## list(mean = mean, median = median)
##
## # Auto named with `tibble::lst()`:
## tibble::lst(mean, median)
##
## # Using lambdas
## list(~ mean(., trim = .2), ~ median(., na.rm = TRUE))## # A tibble: 9 x 7
## popID carinata_mean carinulata_mean sublineata_mean carinata_var
## * <chr> <dbl> <dbl> <dbl> <dbl>
## 1 15TX 0.834 0 0.039 0.00568
## 2 16TX 0.793 0 0.207 NA
## 3 18TX 0.65 0 0.35 NA
## 4 20NM 0.281 0 0.688 0.0797
## 5 21NM 0.137 0 0.748 0.0375
## 6 26NM 0.200 0 0.784 0.0205
## 7 27NM 0.36 0 0.639 NA
## 8 28NM 0.468 0.0999 0.425 0.0397
## 9 6KS 0 0.088 0.912 NA
## # … with 2 more variables: carinulata_var <dbl>, sublineata_var <dbl>Which sites have multiple spp, even without admixture?
foo <- q_popmap %>%
filter_all(any_vars(. > .95)) %>%
group_by(popID) %>%
summarize_at(vars(elongata:carinata), funs(sum)) %>%
merge(sample_size)
max_ancestry <- apply(X=foo[,2:5], MARGIN=1, FUN=max)
### if a cluster has less than foo$sample_size-.5, will return any pop with an ind less than 50% of majority ancestry
mixed_pops <- foo[foo$n - .5 > max_ancestry,]
mixed_pops## popID elongata sublineata carinulata carinata n
## 3 15TX 1.580 0.158 0.005 14.252 16
## 4 16TX 16.996 0.209 0.000 0.794 18
## 5 18TX 0.002 14.249 2.092 0.655 17
## 8 20NM 0.125 15.647 0.003 1.223 17
## 9 21NM 0.328 13.295 0.019 0.355 14
## 10 26NM 0.338 15.372 0.019 2.270 18
## 11 27NM 0.321 19.176 0.002 0.495 20
## 12 28NM 0.094 3.035 8.596 5.274 17
## 32 6KS 0.000 0.915 1.088 7.995 10perc_admixed <- merge(sample_size, num_admixed_inds) %>%
mutate(perc_admixed = num_admixed/n) %>%
arrange(desc(perc_admixed))
perc_admixed## popID n num_admixed perc_admixed
## 1 26NM 18 11 0.61111111
## 2 28NM 17 7 0.41176471
## 3 15TX 16 4 0.25000000
## 4 20NM 17 4 0.23529412
## 5 21NM 14 2 0.14285714
## 6 6KS 10 1 0.10000000
## 7 18TX 17 1 0.05882353
## 8 16TX 18 1 0.05555556
## 9 27NM 20 1 0.05000000mix_admix <- full_join(mixed_pops, perc_admixed[,c('popID', 'perc_admixed')]) %>%
arrange(desc(perc_admixed))## Joining, by = "popID"mix_admix## popID elongata sublineata carinulata carinata n perc_admixed
## 1 26NM 0.338 15.372 0.019 2.270 18 0.61111111
## 2 28NM 0.094 3.035 8.596 5.274 17 0.41176471
## 3 15TX 1.580 0.158 0.005 14.252 16 0.25000000
## 4 20NM 0.125 15.647 0.003 1.223 17 0.23529412
## 5 21NM 0.328 13.295 0.019 0.355 14 0.14285714
## 6 6KS 0.000 0.915 1.088 7.995 10 0.10000000
## 7 18TX 0.002 14.249 2.092 0.655 17 0.05882353
## 8 16TX 16.996 0.209 0.000 0.794 18 0.05555556
## 9 27NM 0.321 19.176 0.002 0.495 20 0.05000000How many 'pure' sites did we detect?
setdiff(popmap$popID, mix_admix$popID )## [1] "11CO" "12AZ" "19TX" "1WY" "2CO"
## [6] "31NV" "32UT" "33NV" "34UT" "37CR"
## [11] "38CR" "39CR" "41CR" "43GR" "44GR"
## [16] "46CH" "47TX" "49TX" "4CO" "50TX"
## [21] "51TX" "52TX" "5CO" "8OK" "CARINA_LAB"
## [26] "SUB_LAB"length(setdiff(popmap$popID, mix_admix$popID ))## [1] 26Pull pure parental species for reference down the road.
pure_parentals <- ci %>%
dplyr::filter(carinata_L > 1 - t_low | elongata_L > 1 - t_low |
sublineata_L > 1 - t_low | carinulata_L > 1 - t_low) %>%
merge(popmap)
pure_parentals$assignment <- colnames(pure_parentals[,c('elongata', 'carinata', 'carinulata', 'sublineata')])[max.col(pure_parentals[,c('elongata', 'carinata', 'carinulata', 'sublineata')])]
table(pure_parentals$assignment)##
## carinata carinulata elongata sublineata
## 64 179 92 166pure_parentals_popmap <- pure_parentals[,c('sampID', 'assignment')]
write.table(pure_parentals_popmap, "../info/pure_parentals_popmap.tsv", col.names = F, quote = F, row.names = F, sep = "\t")Substructure
Subset individuals for substructure of elongata, carinulata, and everything but carinulata ("others").
elongata <- dplyr::filter(ci, elongata_L > t_low) %>%
select('sampID', 'popID')
nrow(elongata)## [1] 99write.table(elongata,
file = "../info/elongata_popmap.tsv",
sep = "\t",
quote=F, col.names = F, row.names = F)
carinulata <- dplyr::filter(ci, carinulata_L > t_low) %>%
select('sampID', 'popID')
nrow(carinulata)## [1] 182write.table(carinulata,
file = "../info/carinulata_popmap.tsv",
sep = "\t",
quote=F, col.names = F, row.names = F)
others_popmap <- anti_join(popmap, carinulata)## Joining, by = c("sampID", "popID")write.table(others_popmap,
file = "../info/others_popmap.tsv",
sep = "\t",
quote=F, col.names = F, row.names = F)Population structure within D. carinulata
sfiles <- list.files(path="../substructure/carinulata/structure_analysis/Results/",
pattern=glob2rx("outfile_*f"),full.names=TRUE)
carinu_slist <- readQ(sfiles,filetype="structure",
indlabfromfile = TRUE,
readci = TRUE)
popmap <- left_join(data.frame(sampID = rownames(carinu_slist[[1]])),
masterpopmap)## Joining, by = "sampID"popID <- data.frame(popID = popmap$popID)
popID <- data.frame(lapply(popID, as.character), stringsAsFactors=FALSE)
sapply(popID, is.character)## popID
## TRUEpop_order <- sample_meta.dat[,c("popID", "Order")]
pop_order <- semi_join(pop_order, popID) %>%
arrange(Order)## Joining, by = "popID"pop_order <- pop_order[,1]
legend_colors## cbbPalette human_colors spp
## 1 #56B4E9 blue carinata
## 2 #E69F00 orange elongata
## 3 #009E73 green sublineata
## 4 #F0E442 yellow carinulata
## 5 #0072B2 dark blue elongata_sub
## 6 #D55E00 red carinulata_sub
## 7 #9370DB purple othercarinu_legend <- legend_colors[c(4,6,7),]
carinu_cbbPalette <- carinu_legend[,'cbbPalette']
carinu_cbbPalette## [1] "#F0E442" "#D55E00" "#9370DB"sr1 <- summariseQ(tabulateQ(carinu_slist))
p <- evannoMethodStructure(data=sr1,exportplot=F,returnplot=T,returndata=F,basesize=12,linesize=0.7)
grid.arrange(p)
evannoMethodStructure(data=sr1,
exportplot=T,
exportpath = "../substructure/carinulata/structure_analysis/plots/",
outputfilename = "carinulata_deltaK",
returndata=F,
height = 10, width = 9)## ../substructure/carinulata/structure_analysis/plots//carinulata_deltaK.png exported.plotQ(alignK(carinu_slist[c(21:30)]),
imgoutput="join",
exportpath = "../substructure/carinulata/structure_analysis/plots/",
outputfilename = "Joined-K2",
exportplot=T,
showyaxis=T,
showsp = F,
clustercol=carinu_cbbPalette,
subsetgrp = pop_order,
height = .1,
ordergrp = T,
grplab = popID,
linepos = 1,
divtype = 1,
grplabangle=45,
grplabpos = .8,
grplabspacer = 0,
grplabheight=70,
width=100,
showtitle=T,
grplabsize = 4,
titlelab="D. carinulata Population Structure, K=2",
titlesize = 20)## Drawing plot ...## ../substructure/carinulata/structure_analysis/plots//Joined-K2.png exported.plotQ(alignK(carinu_slist[c(31:40)]),
imgoutput="join",
exportpath = "../substructure/carinulata/structure_analysis/plots/",
outputfilename = "Joined-K3",
exportplot=T,
showyaxis=T,
showsp = F,
clustercol=carinu_cbbPalette,
subsetgrp = pop_order,
height = .1,
ordergrp = T,
grplab = popID,
linepos = 1,
divtype = 1,
grplabangle=45,
grplabpos = .8,
grplabspacer = 0,
grplabheight=70,
width=100,
showtitle=T,
grplabsize = 4,
titlelab="D. carinulata Population Structure, K=3",
titlesize = 20)## Drawing plot ...## ../substructure/carinulata/structure_analysis/plots//Joined-K3.png exported.k4 <- plotQ(alignK(carinu_slist[c(41:50)]),
imgoutput="join",
returnplot=T,exportplot=F,basesize=11,
clustercol=cbbPalette,
subsetgrp = pop_order,
ordergrp = T,
grplab = popID,
grplabangle=-90)
grid.arrange(k4$plot[[1]])
k5 <- plotQ(alignK(carinu_slist[c(51:60)]),
imgoutput="join",
subsetgrp = pop_order,
returnplot=T,exportplot=F,basesize=11,
clustercol=cbbPalette,
ordergrp = T,
grplab = popID,
grplabangle=-90)
grid.arrange(k5$plot[[1]])
k6 <- plotQ(alignK(carinu_slist[c(61:70)]),
imgoutput="join",
subsetgrp = pop_order,
returnplot=T,exportplot=F,basesize=11,
clustercol=cbbPalette,
ordergrp = T,
grplab = popID,
grplabangle=-90)
grid.arrange(k6$plot[[1]])
k7 <- plotQ(alignK(carinu_slist[c(71:80)]),
imgoutput="join",
subsetgrp = pop_order,
returnplot=T,exportplot=F,basesize=11,
clustercol=cbbPalette,
ordergrp = T,
grplab = popID,
grplabangle=-90)
grid.arrange(k7$plot[[1]])
plotQ(carinu_slist[32],
exportpath = "../substructure/carinulata/structure_analysis/plots/",
outputfilename = "MainTextFigure_carinuK3",
exportplot=T,
showyaxis=T,
clustercol=carinu_cbbPalette,
subsetgrp = pop_order,
sortind="all",
grplab = popID,
linepos = 1,
linesize=0.8,
pointsize = 5,
divtype = 1,
divsize = 1.2,
grplabangle= 90,
grplabpos = .9,
grplabspacer = 0,
grplabheight=40,
grplabsize = 25,
height = .05,
width=100,
showlegend = T,
showsp = F,
legendkeysize=40,legendtextsize=40, legendpos="left", legendlab= c("D. carinulata Fukang Ecotype", "D. carinulata Chilik Ecotype", "Other"),
basesize=40)## Drawing plot ...## ../substructure/carinulata/structure_analysis/plots//MainTextFigure_carinuK3.png exported.Population structure within D. elongata
sfiles <- list.files(path="../substructure/elongata/structure_analysis/Results/",
pattern=glob2rx("outfile_*f"),full.names=TRUE)
elong_slist <- readQ(sfiles,filetype="structure",
indlabfromfile = TRUE,
readci = TRUE)
popmap <- left_join(data.frame(sampID = rownames(elong_slist[[1]])),
masterpopmap)## Joining, by = "sampID"popID <- data.frame(popID = popmap$popID)
popID <- data.frame(lapply(popID, as.character), stringsAsFactors=FALSE)
sapply(popID, is.character)## popID
## TRUEpop_order <- sample_meta.dat[,c("popID", "Order")]
pop_order <- semi_join(pop_order, popID) %>%
arrange(Order)## Joining, by = "popID"pop_order <- pop_order[,1]
legend_colors## cbbPalette human_colors spp
## 1 #56B4E9 blue carinata
## 2 #E69F00 orange elongata
## 3 #009E73 green sublineata
## 4 #F0E442 yellow carinulata
## 5 #0072B2 dark blue elongata_sub
## 6 #D55E00 red carinulata_sub
## 7 #9370DB purple otherelong_legend <- legend_colors[c(2,5,7),]
elong_cbbPalette <- elong_legend[,'cbbPalette']
elong_cbbPalette## [1] "#E69F00" "#0072B2" "#9370DB"sr1 <- summariseQ(tabulateQ(elong_slist))
evannoMethodStructure(data=sr1,
exportplot=T,
exportpath = "../substructure/elongata/structure_analysis/plots/",
outputfilename = "elongata_deltaK",
returndata=F,
height = 10, width = 9)## ../substructure/elongata/structure_analysis/plots//elongata_deltaK.png exported.k2 <- plotQ(alignK(elong_slist[c(21:30)]),
imgoutput="join",
returnplot=T,exportplot=F,basesize=11,
clustercol=cbbPalette,
subsetgrp = pop_order,
ordergrp = T,
grplab = popID,
linepos = 1,
grplabangle=-45,
grplabpos = 0,
grplabheight=20)
grid.arrange(k2$plot[[1]])
k3 <- plotQ(alignK(elong_slist[c(31:40)]),
imgoutput="join",
returnplot=T,exportplot=F,basesize=11,
clustercol=cbbPalette,
ordergrp = T,
grplab = popID,
grplabangle=-90)
grid.arrange(k3$plot[[1]])
k4 <- plotQ(alignK(elong_slist[c(41:50)]),
imgoutput="join",
returnplot=T,exportplot=F,basesize=11,
clustercol=cbbPalette,
ordergrp = T,
grplab = popID,
grplabangle=-90)
grid.arrange(k4$plot[[1]])
k5 <- plotQ(alignK(elong_slist[c(51:60)]),
imgoutput="join",
returnplot=T,exportplot=F,basesize=11,
clustercol=cbbPalette,
ordergrp = T,
grplab = popID,
grplabangle=-90)
grid.arrange(k5$plot[[1]])
k6 <- plotQ(alignK(elong_slist[c(61:70)]),
imgoutput="join",
returnplot=T,exportplot=F,basesize=11,
clustercol=cbbPalette,
ordergrp = T,
grplab = popID,
grplabangle=-90)
grid.arrange(k6$plot[[1]])
k7 <- plotQ(alignK(elong_slist[c(71:80)]),
imgoutput="join",
returnplot=T,exportplot=F,basesize=11,
ordergrp = T,
grplab = popID,
grplabangle=-90)
grid.arrange(k7$plot[[1]])
plotQ(elong_slist[37],
exportpath = "../substructure/elongata/structure_analysis/plots/",
outputfilename = "MainTextFigure_elongataK3",
exportplot=T,
showyaxis=F,
clustercol=elong_cbbPalette,
subsetgrp = pop_order,
sortind="all",
grplab = popID,
linepos = 1,
linesize=0.8,
pointsize = 5,
divtype = 1,
divsize = 1.2,
grplabangle= 90,
grplabpos = .9,
grplabspacer = 0,
grplabheight=40,
grplabsize = 25,
height = .05,
width=100,
showlegend = F,
showsp = F,
basesize=40)## Drawing plot ...## ../substructure/elongata/structure_analysis/plots//MainTextFigure_elongataK3.png exported.Is there a latitudinal gradient found for ecotypes or ancestry?
library(ggtext)
global_q.dat <- slist[[41]]
global_q.dat$sampID <- rownames(global_q.dat)
global_q.dat.popmap.meta <- merge(global_q.dat,
masterpopmap) %>%
merge(sample_meta.dat) %>%
dplyr::filter(Grp == "potential_hybrid_zone")
sub_plot <- ggplot(global_q.dat.popmap.meta, aes(x = Latitude, y = Cluster2)) +
geom_smooth(method = "lm", color = "red") +
geom_jitter(aes(color = popID, shape = cat), alpha = .8, size = 4) +
ylim(c(0,1)) +
ylab("Proportion of *D. sublineata* ancestry") +
guides(col = guide_legend(ncol = 2)) +
theme_bw(base_size = 12)+
theme(axis.title.y = element_markdown())
sub_plot## `geom_smooth()` using formula 'y ~ x'## Warning: Removed 6 rows containing missing values (geom_smooth).## Warning: Removed 83 rows containing missing values (geom_point).
ggsave("sublineata_ancestry_latitude.jpg")## Saving 7 x 5 in image
## `geom_smooth()` using formula 'y ~ x'## Warning: Removed 6 rows containing missing values (geom_smooth).## Warning: Removed 89 rows containing missing values (geom_point).fit_carina <- lm(global_q.dat.popmap.meta$Cluster4 ~ global_q.dat.popmap.meta$Latitude)
summary(fit_carina)##
## Call:
## lm(formula = global_q.dat.popmap.meta$Cluster4 ~ global_q.dat.popmap.meta$Latitude)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.81015 -0.26457 0.06703 0.08528 0.67864
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -3.042344 0.251300 -12.11 <2e-16 ***
## global_q.dat.popmap.meta$Latitude 0.101408 0.007729 13.12 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3187 on 275 degrees of freedom
## Multiple R-squared: 0.385, Adjusted R-squared: 0.3827
## F-statistic: 172.1 on 1 and 275 DF, p-value: < 2.2e-16fit_sub <- lm(global_q.dat.popmap.meta$Cluster2 ~ global_q.dat.popmap.meta$Latitude)
summary(fit_sub)##
## Call:
## lm(formula = global_q.dat.popmap.meta$Cluster2 ~ global_q.dat.popmap.meta$Latitude)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.75212 -0.12038 -0.01282 0.24788 0.94804
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.570464 0.270919 16.87 <2e-16 ***
## global_q.dat.popmap.meta$Latitude -0.121256 0.008333 -14.55 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3436 on 275 degrees of freedom
## Multiple R-squared: 0.435, Adjusted R-squared: 0.433
## F-statistic: 211.8 on 1 and 275 DF, p-value: < 2.2e-16plot(fit_sub, which=1, col=c("blue"))
plot(fit_sub, which=2, col=c("blue"))
global_q.dat.popmap.meta$predicted <- predict(fit_sub) # Save the predicted values
global_q.dat.popmap.meta$residuals <- residuals(fit_sub) # Save the residual values
ggplot(global_q.dat.popmap.meta, aes(x = Latitude, y = Cluster2)) +
geom_smooth(method = "lm", se = FALSE, color = "lightgrey") + # regression line
geom_segment(aes(xend = Latitude, yend = predicted), alpha = .2) + # draw line from point to line
geom_point(aes(color = abs(residuals), size = abs(residuals))) + # size of the points
scale_color_continuous(low = "green", high = "red") + # colour of the points mapped to residual size - green smaller, red larger
guides(color = FALSE, size = FALSE) + # Size legend removed
geom_point(aes(y = predicted), shape = 1) +
ylab("Proportion of D. sublineata ancestry") +
theme_bw()## `geom_smooth()` using formula 'y ~ x'
ggsave("sublineata_ancestry_latitude_fit.jpg")## Saving 7 x 5 in image
## `geom_smooth()` using formula 'y ~ x'carinu_q.dat <- carinu_slist[[32]]
carinu_q.dat$sampID <- rownames(carinu_q.dat)
carinu_q.popmap.meta <- merge(carinu_q.dat,
masterpopmap) %>%
merge(sample_meta.dat) #%>%
# dplyr::filter(Cluster3 < .5) # %>%
# dplyr::filter( popID != "46CH" & popID != "1WY"& popID != "18TX" & popID != "28NM")
carinu_plot <- ggplot(carinu_q.popmap.meta) +
geom_smooth(aes(Latitude, Cluster2), method = "lm", color = "red") +
geom_jitter(aes(x = Latitude, y = Cluster2, color = popID, shape = cat), alpha = .8, size = 4) +
ylab("Proportion of Chilik ancestry") +
ylim(c(0,1)) +
guides(col = guide_legend(ncol = 2)) +
theme_bw(base_size = 12)
carinu_plot## `geom_smooth()` using formula 'y ~ x'## Warning: Removed 2 rows containing missing values (geom_smooth).## Warning: Removed 1 rows containing missing values (geom_point).
ggsave("fukang_ancestry_latitude.jpg")## Saving 7 x 5 in image
## `geom_smooth()` using formula 'y ~ x'## Warning: Removed 2 rows containing missing values (geom_smooth).fit_carinu <- lm(carinu_q.popmap.meta$Cluster2 ~ carinu_q.popmap.meta$Latitude)
summary(fit_carinu)##
## Call:
## lm(formula = carinu_q.popmap.meta$Cluster2 ~ carinu_q.popmap.meta$Latitude)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.86712 -0.47929 0.02838 0.40506 0.56583
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.641788 0.297754 8.872 7.03e-16 ***
## carinu_q.popmap.meta$Latitude -0.056325 0.007543 -7.468 3.39e-12 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.4096 on 180 degrees of freedom
## Multiple R-squared: 0.2365, Adjusted R-squared: 0.2323
## F-statistic: 55.77 on 1 and 180 DF, p-value: 3.394e-12plot(fit_carinu, which=1, col=c("blue"))
plot(fit_carinu, which=2, col=c("blue"))
carinu_q.popmap.meta$predicted <- predict(fit_carinu) # Save the predicted values
carinu_q.popmap.meta$residuals <- residuals(fit_carinu) # Save the residual values
ggplot(carinu_q.popmap.meta, aes(x = Latitude, y = Cluster2)) +
geom_smooth(method = "lm", se = FALSE, color = "lightgrey") + # regression line
geom_segment(aes(xend = Latitude, yend = predicted), alpha = .2) + # draw line from point to line
geom_point(aes(color = abs(residuals), size = abs(residuals))) + # size of the points
scale_color_continuous(low = "green", high = "red") +
ylab("Proportion of Fukang ancestry") + # colour of the points mapped to residual size - green smaller, red larger
guides(color = FALSE, size = FALSE) +
ylab("Proportion of Chilik ancestry") +# Size legend removed
geom_point(aes(y = predicted), shape = 1) +
theme_bw()## `geom_smooth()` using formula 'y ~ x'
ggsave("fukang_ancestry_latitude_fit.jpg")## Saving 7 x 5 in image
## `geom_smooth()` using formula 'y ~ x'library(cowplot)
plots <- align_plots( sub_plot, carinu_plot, align = 'v', axis = 'l')## `geom_smooth()` using formula 'y ~ x'## Warning: Removed 6 rows containing missing values (geom_smooth).## Warning: Removed 86 rows containing missing values (geom_point).## `geom_smooth()` using formula 'y ~ x'## Warning: Removed 2 rows containing missing values (geom_smooth).plot_grid(plots[[1]],
plots[[2]],
ncol = 1,
labels = c('A', 'B'))
ggsave("sublineata_carinulata_ancestry_latitude.png", dpi = 300, height = 8, width = 8, units = "in")